home *** CD-ROM | disk | FTP | other *** search
/ Aminet 5 / Aminet 5 - March 1995.iso / Aminet / game / role / Ang261Lib.lha / src / random.c < prev    next >
C/C++ Source or Header  |  1994-10-22  |  13KB  |  375 lines
  1. /*
  2.  * Copyright (c) 1983 Regents of the University of California. All rights
  3.  * reserved. 
  4.  *
  5.  * Redistribution and use in source and binary forms are permitted provided that
  6.  * the above copyright notice and this paragraph are duplicated in all such
  7.  * forms and that any documentation, advertising materials, and other
  8.  * materials related to such distribution and use acknowledge that the
  9.  * software was developed by the University of California, Berkeley.  The
  10.  * name of the University may not be used to endorse or promote products
  11.  * derived from this software without specific prior written permission. THIS
  12.  * SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR IMPLIED
  13.  * WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED WARRANTIES OF
  14.  * MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE. 
  15.  */
  16.  
  17. /*
  18.  * Some minor changes were made by me (C. Swiger): 
  19.  *
  20.  * removed SCCS gunk added external prototypes for both ANSI & non-ANSI
  21.  * compilers removed unused variable "j" from srandom() re-indented and
  22.  * respaced the code -- someone has one ugly coding style wrt whitespace
  23.  * ;-).... 
  24.  */
  25.  
  26. #include <stdio.h>
  27.  
  28. /* some prototypes - cba */
  29. #ifndef NO_LINT_ARGS
  30. #ifdef __STDC__
  31. extern long random(void);
  32. extern void srandom(unsigned int);
  33. extern char *initstate(unsigned int, char *, int);
  34. extern char *setstate(char *);
  35. #else
  36. extern long random();
  37. extern int srandom();
  38. extern char *initstate();
  39. extern char *setstate();
  40. #endif
  41. #endif
  42.  
  43. /*
  44.  * random.c: An improved random number generation package.  In addition to
  45.  * the standard rand()/srand() like interface, this package also has a
  46.  * special state info interface.  The initstate() routine is called with a
  47.  * seed, an array of bytes, and a count of how many bytes are being passed
  48.  * in; this array is then initialized to contain information for random
  49.  * number generation with that much state information.  Good sizes for the
  50.  * amount of state information are 32, 64, 128, and 256 bytes.  The state can
  51.  * be switched by calling the setstate() routine with the same array as was
  52.  * initiallized with initstate(). By default, the package runs with 128 bytes
  53.  * of state information and generates far better random numbers than a linear
  54.  * congruential generator. If the amount of state information is less than 32
  55.  * bytes, a simple linear congruential R.N.G. is used. Internally, the state
  56.  * information is treated as an array of longs; the zeroeth element of the
  57.  * array is the type of R.N.G. being used (small integer); the remainder of
  58.  * the array is the state information for the R.N.G.  Thus, 32 bytes of state
  59.  * information will give 7 longs worth of state information, which will allow
  60.  * a degree seven polynomial.  (Note: the zeroeth word of state information
  61.  * also has some other information stored in it -- see setstate() for
  62.  * details). The random number generation technique is a linear feedback
  63.  * shift register approach, employing trinomials (since there are fewer terms
  64.  * to sum up that way).  In this approach, the least significant bit of all
  65.  * the numbers in the state table will act as a linear feedback shift
  66.  * register, and will have period 2^deg - 1 (where deg is the degree of the
  67.  * polynomial being used, assuming that the polynomial is irreducible and
  68.  * primitive).  The higher order bits will have longer periods, since their
  69.  * values are also influenced by pseudo-random carries out of the lower bits. 
  70.  * The total period of the generator is approximately deg*(2**deg - 1); thus
  71.  * doubling the amount of state information has a vast influence on the
  72.  * period of the generator. Note: the deg*(2**deg - 1) is an approximation
  73.  * only good for large deg, when the period of the shift register is the
  74.  * dominant factor.  With deg equal to seven, the period is actually much
  75.  * longer than the 7*(2**7 - 1) predicted by this formula. 
  76.  */
  77.  
  78.  
  79. /*
  80.  * For each of the currently supported random number generators, we have a
  81.  * break value on the amount of state information (you need at least this
  82.  * many bytes of state info to support this random number generator), a
  83.  * degree for the polynomial (actually a trinomial) that the R.N.G. is based
  84.  * on, and the separation between the two lower order coefficients of the
  85.  * trinomial. 
  86.  */
  87.  
  88. #define TYPE_0          0  /* linear congruential */
  89. #define BREAK_0         8
  90. #define DEG_0           0
  91. #define SEP_0           0
  92.  
  93. #define TYPE_1          1  /* x**7 + x**3 + 1 */
  94. #define BREAK_1         32
  95. #define DEG_1           7
  96. #define SEP_1           3
  97.  
  98. #define TYPE_2          2  /* x**15 + x + 1 */
  99. #define BREAK_2         64
  100. #define DEG_2           15
  101. #define SEP_2           1
  102.  
  103. #define TYPE_3          3  /* x**31 + x**3 + 1 */
  104. #define BREAK_3         128
  105. #define DEG_3           31
  106. #define SEP_3           3
  107.  
  108. #define TYPE_4          4  /* x**63 + x + 1 */
  109. #define BREAK_4         256
  110. #define DEG_4           63
  111. #define SEP_4           1
  112.  
  113.  
  114. /*
  115.  * Array versions of the above information to make code run faster -- relies
  116.  * on fact that TYPE_i == i. 
  117.  */
  118.  
  119. #define MAX_TYPES       5  /* max number of types above */
  120.  
  121. static int degrees[MAX_TYPES] = {DEG_0, DEG_1, DEG_2, DEG_3, DEG_4};
  122.  
  123. static int seps[MAX_TYPES] = {SEP_0, SEP_1, SEP_2, SEP_3, SEP_4};
  124.  
  125.  
  126.  
  127. /*
  128.  * Initially, everything is set up as if from : initstate( 1, &randtbl, 128);
  129.  * Note that this initialization takes advantage of the fact that
  130.  * srandom() advances the front and rear pointers 10*rand_deg times, and
  131.  * hence the rear pointer which starts at 0 will also end up at zero; thus
  132.  * the zeroeth element of the state information, which contains info about
  133.  * the current position of the rear pointer is just MAX_TYPES*(rptr - state)
  134.  * + TYPE_3 == TYPE_3. 
  135.  */
  136.  
  137. static long         randtbl[DEG_3 + 1] = {TYPE_3,
  138.                  0x9a319039, 0x32d9c024, 0x9b663182, 0x5da1f342,
  139.                  0xde3b81e0, 0xdf0a6fb5, 0xf103bc02, 0x48f340fb,
  140.                  0x7449e56b, 0xbeb1dbb0, 0xab5c5918, 0x946554fd,
  141.                  0x8c2e680f, 0xeb3d799f, 0xb11ee0b7, 0x2d436b86,
  142.                  0xda672e2a, 0x1588ca88, 0xe369735d, 0x904f35f7,
  143.                  0xd7158fd6, 0x6fa6f051, 0x616e6b96, 0xac94efdc,
  144.                  0x36413f93, 0xc622c298, 0xf5a42ab8, 0x8a88d77b,
  145.                     0xf5ad9d0e, 0x8999220b, 0x27fb47b9};
  146.  
  147. /*
  148.  * fptr and rptr are two pointers into the state info, a front and a rear
  149.  * pointer.  These two pointers are always rand_sep places aparts, as they
  150.  * cycle cyclically through the state information.  (Yes, this does mean we
  151.  * could get away with just one pointer, but the code for random() is more
  152.  * efficient this way).  The pointers are left positioned as they would be
  153.  * from the call initstate( 1, randtbl, 128 ) (The position of the rear
  154.  * pointer, rptr, is really 0 (as explained above in the initialization of
  155.  * randtbl) because the state table pointer is set to point to randtbl[1] (as
  156.  * explained below). 
  157.  */
  158.  
  159. static long        *fptr = &randtbl[SEP_3 + 1];
  160. static long        *rptr = &randtbl[1];
  161.  
  162.  
  163.  
  164. /*
  165.  * The following things are the pointer to the state information table, the
  166.  * type of the current generator, the degree of the current polynomial being
  167.  * used, and the separation between the two pointers. Note that for
  168.  * efficiency of random(), we remember the first location of the state
  169.  * information, not the zeroeth.  Hence it is valid to access state[-1],
  170.  * which is used to store the type of the R.N.G. Also, we remember the last
  171.  * location, since this is more efficient than indexing every time to find
  172.  * the address of the last element to see if the front and rear pointers have
  173.  * wrapped. 
  174.  */
  175.  
  176. static long        *state = &randtbl[1];
  177.  
  178. static int          rand_type = TYPE_3;
  179. static int          rand_deg = DEG_3;
  180. static int          rand_sep = SEP_3;
  181.  
  182. static long        *end_ptr = &randtbl[DEG_3 + 1];
  183.  
  184.  
  185.  
  186. /*
  187.  * srandom: Initialize the random number generator based on the given seed. 
  188.  * If the type is the trivial no-state-information type, just remember the
  189.  * seed. Otherwise, initializes state[] based on the given "seed" via a
  190.  * linear congruential generator.  Then, the pointers are set to known
  191.  * locations that are exactly rand_sep places apart.  Lastly, it cycles the
  192.  * state information a given number of times to get rid of any initial
  193.  * dependencies introduced by the L.C.R.N.G. Note that the initialization of
  194.  * randtbl[] for default usage relies on values produced by this routine. 
  195.  */
  196.  
  197. #ifdef __STDC__
  198. void
  199. #else
  200. int
  201. #endif /* __STDC__ */
  202. srandom(x)
  203.     unsigned int        x;
  204. {
  205.     register int        i;
  206.  
  207.     if (rand_type == TYPE_0) {
  208.     state[0] = x;
  209.     } else {
  210.     state[0] = x;
  211.     for (i = 1; i < rand_deg; i++) {
  212.         state[i] = 1103515245 * state[i - 1] + 12345;
  213.     }
  214.     fptr = &state[rand_sep];
  215.     rptr = &state[0];
  216.     for (i = 0; i < 10 * rand_deg; i++)
  217.         random();
  218.     }
  219. }
  220.  
  221.  
  222. /*
  223.  * initstate: Initialize the state information in the given array of n bytes
  224.  * for future random number generation.  Based on the number of bytes we are
  225.  * given, and the break values for the different R.N.G.'s, we choose the best
  226.  * (largest) one we can and set things up for it.  srandom() is then called
  227.  * to initialize the state information. Note that on return from srandom(),
  228.  * we set state[-1] to be the type multiplexed with the current value of the
  229.  * rear pointer; this is so successive calls to initstate() won't lose this
  230.  * information and will be able to restart with setstate(). Note: the first
  231.  * thing we do is save the current state, if any, just like setstate() so
  232.  * that it doesn't matter when initstate is called. Returns a pointer to the
  233.  * old state. 
  234.  */
  235.  
  236. char               *
  237. initstate(seed, arg_state, n)
  238.     unsigned            seed;       /* seed for R. N. G. */
  239.     char               *arg_state; /* pointer to state array */
  240.     int                 n;       /* # bytes of state info */
  241. {
  242.     register char      *ostate = (char *)(&state[-1]);
  243.  
  244.     if (rand_type == TYPE_0)
  245.     state[-1] = rand_type;
  246.     else
  247.     state[-1] = MAX_TYPES * (rptr - state) + rand_type;
  248.     if (n < BREAK_1) {
  249.     if (n < BREAK_0) {
  250.         fprintf(stderr,
  251.             "initstate: not enough state (%d bytes) with which to do jack; ignored.", n);
  252.         return (char *)0;
  253.     }
  254.     rand_type = TYPE_0;
  255.     rand_deg = DEG_0;
  256.     rand_sep = SEP_0;
  257.     } else {
  258.     if (n < BREAK_2) {
  259.         rand_type = TYPE_1;
  260.         rand_deg = DEG_1;
  261.         rand_sep = SEP_1;
  262.     } else {
  263.         if (n < BREAK_3) {
  264.         rand_type = TYPE_2;
  265.         rand_deg = DEG_2;
  266.         rand_sep = SEP_2;
  267.         } else {
  268.         if (n < BREAK_4) {
  269.             rand_type = TYPE_3;
  270.             rand_deg = DEG_3;
  271.             rand_sep = SEP_3;
  272.         } else {
  273.             rand_type = TYPE_4;
  274.             rand_deg = DEG_4;
  275.             rand_sep = SEP_4;
  276.         }
  277.         }
  278.     }
  279.     }
  280.     state = &(((long *)arg_state)[1]);    /* first location */
  281.     end_ptr = &state[rand_deg];       /* must set end_ptr before srandom */
  282.     srandom(seed);
  283.     if (rand_type == TYPE_0)
  284.     state[-1] = rand_type;
  285.     else
  286.     state[-1] = MAX_TYPES * (rptr - state) + rand_type;
  287.     return (ostate);
  288. }
  289.  
  290.  
  291.  
  292. /*
  293.  * setstate: Restore the state from the given state array. Note: it is
  294.  * important that we also remember the locations of the pointers in the
  295.  * current state information, and restore the locations of the pointers from
  296.  * the old state information.  This is done by multiplexing the pointer
  297.  * location into the zeroeth word of the state information. Note that due to
  298.  * the order in which things are done, it is OK to call setstate() with the
  299.  * same state as the current state. Returns a pointer to the old state
  300.  * information. 
  301.  */
  302.  
  303. char               *
  304. setstate(arg_state)
  305.     char               *arg_state;
  306. {
  307.     register long      *new_state = (long *)arg_state;
  308.     register int        type = new_state[0] % MAX_TYPES;
  309.     register int        rear = new_state[0] / MAX_TYPES;
  310.     char               *ostate = (char *)(&state[-1]);
  311.  
  312.     if (rand_type == TYPE_0)
  313.     state[-1] = rand_type;
  314.     else
  315.     state[-1] = MAX_TYPES * (rptr - state) + rand_type;
  316.     switch (type) {
  317.       case TYPE_0:
  318.       case TYPE_1:
  319.       case TYPE_2:
  320.       case TYPE_3:
  321.       case TYPE_4:{
  322.         rand_type = type;
  323.         rand_deg = degrees[type];
  324.         rand_sep = seps[type];
  325.         break;
  326.     }
  327.       default:
  328.     fprintf(stderr, "setstate: state info has been munged (%d); not changed.", type);
  329.     break;
  330.     }
  331.     state = &new_state[1];
  332.     if (rand_type != TYPE_0) {
  333.     rptr = &state[rear];
  334.     fptr = &state[(rear + rand_sep) % rand_deg];
  335.     }
  336.     end_ptr = &state[rand_deg];       /* set end_ptr too */
  337.     return (ostate);
  338. }
  339.  
  340.  
  341.  
  342. /*
  343.  * random: If we are using the trivial TYPE_0 R.N.G., just do the old linear
  344.  * congruential bit.  Otherwise, we do our fancy trinomial stuff, which is
  345.  * the same in all ther other cases due to all the global variables that have
  346.  * been set up.  The basic operation is to add the number at the rear pointer
  347.  * into the one at the front pointer.  Then both pointers are advanced to the
  348.  * next location cyclically in the table.  The value returned is the sum
  349.  * generated, reduced to 31 bits by throwing away the "least random" low bit.
  350.  * Note: the code takes advantage of the fact that both the front and rear
  351.  * pointers can't wrap on the same call by not testing the rear pointer if
  352.  * the front one has wrapped. Returns a 31-bit random number. 
  353.  */
  354.  
  355. long
  356. random()
  357. {
  358.     long                i;
  359.  
  360.     if (rand_type == TYPE_0) {
  361.     i = state[0] = (state[0] * 1103515245 + 12345) & 0x7fffffff;
  362.     } else {
  363.     *fptr += *rptr;
  364.     i = (*fptr >> 1) & 0x7fffffff;    /* chucking least random bit */
  365.     if (++fptr >= end_ptr) {
  366.         fptr = state;
  367.         ++rptr;
  368.     } else {
  369.         if (++rptr >= end_ptr)
  370.         rptr = state;
  371.     }
  372.     }
  373.     return (i);
  374. }
  375.